home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
Turnbull China Bikeride
/
Turnbull China Bikeride - Disc 1.iso
/
ARGONET
/
PD
/
MATHS
/
RLAB
/
RLAB125.ZIP
/
!RLaB
/
help_ai
/
filter
< prev
next >
Wrap
Text File
|
1994-11-06
|
3KB
|
111 lines
Synopsis: Digital filter implementation
Syntax: filter ( B , A , X )
filter ( B , A , X , Zi )
Description:
Filter is an implementation of the standard difference
equation:
y[n] = b(1)*x[n] + b(2)*x[n-1] + ... b(nb+1)*x[n-nb]
- a(2)*y[n-1] - ... a(na+1)*y[n-na]
The filter is implemented using a method described as a
"Direct Form II Transposed" filter. More for information see
Chapter 6 of "Discrete-Time Signal Processing" by Oppenheim
and Schafer.
The inputs to filter are:
B The numerator coefficients, or zeros of the
system transfer function. The coefficients are
specified in a vector like:
[ b(1) , b(2) , ... b(nb) ]
A The denominator coefficients, or the poles of
the system transfer function. the coefficients
are specified in a vector like:
[ a(1) , a(2) , ... a(na) ]
X A vector of the filter inputs.
Zi [Optional] The initial delays of the filter.
The filter outputs are in a list with element names:
y The filter output. y is a vector of the same
dimension as X.
zf A vector of the final values of the filter
delays.
The A(1) coefficient must be non-zero, as the other
coefficients are divided by A(1).
Below is an implementation of filter() in a r-file - it is
provided for informational usage only.
#--------------------------------------------------------------
# Simplistic version of RLaB's builtin function filter()
# Y = filter ( b, a, x )
# Y = filter ( b, a, x, zi )
#
rfilter = function ( b , a , x , zi )
{
local ( b , a , x , zi )
ntotal = x.nr * x.nc;
M = b.nr * b.nc;
N = a.nr * a.nc;
NN = max ([ M, N ]);
y = zeros (x.nr, x.nc);
# Fix up pole and zero vectors.
# Make them the same length, this makes
# filter's job much easier.
if (N < NN) { a[NN] = 0; }
if (M < NN) { b[NN] = 0; }
# Adjust filter coefficients
if (a[1] == 0) { error ("rfilter: 1st A term must be non-zero"); }
a[2:NN] = a[2:NN] ./ a[1];
b = b ./ a[1];
# Create delay vectors and load inital delays.
# Add an extra term to vi[] to make filter's
# job a little easier. This extra term will
# always be zero.
v = zeros (NN, 1);
vi = zeros (NN+1, 1);
if (exist (zi))
{
vi[1:NN] = zi;
}
#
# Do the work...
#
for (n in 1:ntotal)
{
v[1] = b[1]*x[n] + vi[2];
y[n] = v[1];
for (k in 2:NN)
{
v[k] = b[k]*x[n] - a[k]*v[1] + vi[k+1];
vi[k] = v[k];
}
}
return << y = y; zf = v >>;
};
#--------------------------------------------------------------